{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import geopandas as gpd\n",
    "import pandas as pd\n",
    "import networkx as nx\n",
    "import osmnx as ox\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# 设置代理地址和端口\n",
    "proxy_address = 'http://127.0.0.1:10890'\n",
    "\n",
    "# 设置环境变量\n",
    "os.environ['HTTP_PROXY'] = proxy_address\n",
    "os.environ['HTTPS_PROXY'] = proxy_address\n",
    "\n",
    "# 忽略警告信息\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')\n",
    "\n",
    "# 设置绘图样式\n",
    "%matplotlib inline\n",
    "plt.style.use('ggplot')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 读取Shapefile文件\n",
    "gdf = gpd.read_file('./test_shp/hankou.shp')\n",
    "\n",
    "# 确保坐标参考系为WGS84\n",
    "gdf = gdf.to_crs(epsg=4326)\n",
    "\n",
    "# 提取“Residential”多边形\n",
    "residential_gdf = gdf[gdf['euluc_en'] == 'Residential'].copy()\n",
    "\n",
    "# 定义目标类型列表（排除\"Residential\"和\"Industrial\"）\n",
    "target_types = [\n",
    "    \"Business office\",\n",
    "    \"Commercial service\",\n",
    "    \"Transportation stations\",\n",
    "    \"Administrative\",\n",
    "    \"Educational\",\n",
    "    \"Medical\",\n",
    "    \"Sport and cultural\",\n",
    "    \"Park and greenspace\"\n",
    "]\n",
    "\n",
    "# 提取目标多边形，按类型分类\n",
    "target_gdfs = {}\n",
    "for t_type in target_types:\n",
    "    target_gdf = gdf[gdf['euluc_en'] == t_type].copy()\n",
    "    if not target_gdf.empty:\n",
    "        target_gdfs[t_type] = target_gdf\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "from matplotlib import rcParams\n",
    "\n",
    "# 设置字体为Arial\n",
    "rcParams['font.family'] = 'Arial'\n",
    "\n",
    "# 创建一个绘图对象\n",
    "fig, ax = plt.subplots(figsize=(12, 12))\n",
    "\n",
    "# 绘制所有多边形\n",
    "gdf.plot(ax=ax, color='lightgrey', edgecolor='white')\n",
    "\n",
    "# 绘制“Residential”多边形\n",
    "residential_gdf.plot(ax=ax, color='#fcecb6', edgecolor='white', label='Residential')\n",
    "\n",
    "# 使用不同的颜色绘制各个目标类型的多边形\n",
    "color_map = {\n",
    "    'Residential': '#fcecb6',\n",
    "    'Business office': '#f4c1bf',\n",
    "    'Commercial service': '#eb8682',\n",
    "    'Transportation stations': '#eaebdc',\n",
    "    'Administrative': '#96d5c4',\n",
    "    'Educational': '#bdd2ff',\n",
    "    'Medical': '#a4d3c5',\n",
    "    'Sport and cultural': '#87b0f2',\n",
    "    'Park and greenspace': '#cdedab',\n",
    "}\n",
    "\n",
    "for idx, (t_type, target_gdf) in enumerate(target_gdfs.items()):\n",
    "    target_gdf.plot(ax=ax, color=color_map[t_type], edgecolor='white', label=t_type)\n",
    "\n",
    "# 添加图例\n",
    "ax.legend()\n",
    "# ax.set_title('Residential and Target Areas')\n",
    "ax.axis('off')\n",
    "\n",
    "# 保存为PDF格式，确保矢量化输出\n",
    "plt.savefig('landuse.pdf', format='pdf', bbox_inches='tight')\n",
    "\n",
    "# 显示图形\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 计算“Residential”质心\n",
    "residential_gdf['centroid'] = residential_gdf.geometry.centroid\n",
    "\n",
    "# 计算目标类型质心\n",
    "for t_type, target_gdf in target_gdfs.items():\n",
    "    target_gdf['centroid'] = target_gdf.geometry.centroid\n",
    "    target_gdfs[t_type] = target_gdf\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "from matplotlib import rcParams\n",
    "\n",
    "# 设置字体为Arial\n",
    "rcParams['font.family'] = 'Arial'\n",
    "\n",
    "# 创建一个绘图对象\n",
    "fig, ax = plt.subplots(figsize=(12, 12))\n",
    "\n",
    "# 绘制“Residential”多边形\n",
    "residential_gdf.plot(ax=ax, color='#fcecb6', edgecolor='white', label='Residential', alpha=0.75)\n",
    "\n",
    "# 绘制“Residential”质心\n",
    "residential_gdf.centroid.plot(ax=ax, color='red', markersize=20, label='Residential Centroids')\n",
    "\n",
    "# 绘制目标多边形和质心\n",
    "for idx, (t_type, target_gdf) in enumerate(target_gdfs.items()):\n",
    "    target_gdf.plot(ax=ax, color=color_map[t_type], edgecolor='white', label=t_type, alpha=0.75)\n",
    "    target_gdf.centroid.plot(ax=ax, color='blue', markersize=20, label='Other Centroids')\n",
    "\n",
    "# 获取当前图例的句柄和标签\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "\n",
    "# 去重，保持顺序不变\n",
    "unique_labels = dict(zip(labels, handles))\n",
    "\n",
    "# 添加去重后的图例\n",
    "ax.legend(unique_labels.values(), unique_labels.keys())\n",
    "\n",
    "# 添加标题和关闭坐标轴\n",
    "ax.set_title('Centroids of Residential and Target Areas')\n",
    "ax.axis('off')\n",
    "\n",
    "# 保存为PDF格式，确保矢量化输出\n",
    "plt.savefig('centroid_plot.pdf', format='pdf', bbox_inches='tight')\n",
    "\n",
    "# 显示图形\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import osmnx as ox\n",
    "from matplotlib import rcParams\n",
    "\n",
    "# 设置字体为Arial\n",
    "rcParams['font.family'] = 'Arial'\n",
    "\n",
    "# 计算研究区域的边界\n",
    "minx, miny, maxx, maxy = gdf.total_bounds\n",
    "\n",
    "# 添加缓冲区以覆盖边界外的区域\n",
    "buffer = 0.01  # 可以调整缓冲区的大小\n",
    "north, south, east, west = maxy + buffer, miny - buffer, maxx + buffer, minx - buffer\n",
    "\n",
    "# 获取适用于步行的道路网络\n",
    "G = ox.graph_from_bbox(north, south, east, west, network_type='walk')\n",
    "\n",
    "# 创建绘图对象\n",
    "fig, ax = plt.subplots(figsize=(12, 12))\n",
    "\n",
    "# 绘制道路网络\n",
    "ox.plot_graph(G, ax=ax, node_size=0, edge_color='black', show=False, close=False)\n",
    "\n",
    "# 绘制“Residential”多边形\n",
    "residential_gdf.plot(ax=ax, color='#fcecb6', alpha=0.75, edgecolor='white', label='Residential')\n",
    "\n",
    "# 绘制目标多边形\n",
    "for idx, (t_type, target_gdf) in enumerate(target_gdfs.items()):\n",
    "    target_gdf.plot(ax=ax, color=color_map[t_type], alpha=0.75, edgecolor='white', label=t_type)\n",
    "\n",
    "# 获取当前图例的句柄和标签\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "\n",
    "# 去重，保持顺序不变\n",
    "unique_labels = dict(zip(labels, handles))\n",
    "\n",
    "# 添加去重后的图例\n",
    "ax.legend(unique_labels.values(), unique_labels.keys())\n",
    "\n",
    "# 添加标题并关闭坐标轴\n",
    "ax.set_title('Road Network with Residential and Target Areas')\n",
    "ax.axis('off')\n",
    "\n",
    "# 保存为PDF格式，确保矢量化输出\n",
    "plt.savefig('road_network_plot.pdf', format='pdf', bbox_inches='tight')\n",
    "\n",
    "# 显示图形\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 为每条边添加长度（如果尚未添加）\n",
    "G = ox.distance.add_edge_lengths(G)\n",
    "\n",
    "# 设置平均步行速度（米/秒），例如4.8 km/h\n",
    "walking_speed_mps = 4.8 / 3.6  # 转换为米/秒\n",
    "\n",
    "# 计算每条边的旅行时间\n",
    "for u, v, data in G.edges(data=True):\n",
    "    data['travel_time'] = data['length'] / walking_speed_mps"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import osmnx as ox\n",
    "from matplotlib import rcParams\n",
    "\n",
    "# 设置字体为Arial\n",
    "rcParams['font.family'] = 'Arial'\n",
    "\n",
    "# 为“Residential”质心找到最近的节点\n",
    "residential_gdf['node'] = residential_gdf['centroid'].apply(\n",
    "    lambda point: ox.distance.nearest_nodes(G, X=point.x, Y=point.y)\n",
    ")\n",
    "\n",
    "# 为每个目标类型的质心找到最近的节点\n",
    "for t_type, target_gdf in target_gdfs.items():\n",
    "    target_gdf['node'] = target_gdf['centroid'].apply(\n",
    "        lambda point: ox.distance.nearest_nodes(G, X=point.x, Y=point.y)\n",
    "    )\n",
    "    target_gdfs[t_type] = target_gdf\n",
    "\n",
    "# 获取节点的坐标\n",
    "residential_nodes = G.nodes\n",
    "residential_gdf['node_x'] = residential_gdf['node'].apply(lambda n: residential_nodes[n]['x'])\n",
    "residential_gdf['node_y'] = residential_gdf['node'].apply(lambda n: residential_nodes[n]['y'])\n",
    "\n",
    "for t_type, target_gdf in target_gdfs.items():\n",
    "    target_nodes = G.nodes\n",
    "    target_gdf['node_x'] = target_gdf['node'].apply(lambda n: target_nodes[n]['x'])\n",
    "    target_gdf['node_y'] = target_gdf['node'].apply(lambda n: target_nodes[n]['y'])\n",
    "    target_gdfs[t_type] = target_gdf\n",
    "\n",
    "# 创建绘图对象\n",
    "fig, ax = plt.subplots(figsize=(12, 12))\n",
    "\n",
    "# 绘制“Residential”多边形\n",
    "residential_gdf.plot(ax=ax, color='#fcecb6', alpha=0.75, edgecolor='white', label='Residential')\n",
    "\n",
    "# 绘制道路网络\n",
    "ox.plot_graph(G, ax=ax, node_size=0, edge_color='lightgrey', show=False, close=False)\n",
    "\n",
    "# 绘制“Residential”质心和对应的节点\n",
    "ax.scatter(residential_gdf['centroid'].x, residential_gdf['centroid'].y, color='red', label='Residential Centroids', s=50)\n",
    "ax.scatter(residential_gdf['node_x'], residential_gdf['node_y'], color='#b0d345', label='Residential Nodes', s=30)\n",
    "\n",
    "# 绘制连接线\n",
    "for idx, row in residential_gdf.iterrows():\n",
    "    x_values = [row['centroid'].x, row['node_x']]\n",
    "    y_values = [row['centroid'].y, row['node_y']]\n",
    "    ax.plot(x_values, y_values, color='black', linewidth=1, alpha=0.5)\n",
    "\n",
    "# 获取当前图例的句柄和标签\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "\n",
    "# 去重，保持顺序不变\n",
    "unique_labels = dict(zip(labels, handles))\n",
    "\n",
    "# 添加去重后的图例\n",
    "ax.legend(unique_labels.values(), unique_labels.keys())\n",
    "\n",
    "# 添加标题并关闭坐标轴\n",
    "ax.set_title('Mapping of Residential Centroids to Network Nodes')\n",
    "ax.axis('off')\n",
    "\n",
    "# 保存为PDF格式，确保矢量化输出\n",
    "plt.savefig('centroid_node_mapping.pdf', format='pdf', bbox_inches='tight')\n",
    "\n",
    "# 显示图形\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "res_id = 130"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import networkx as nx\n",
    "import osmnx as ox\n",
    "import geopandas as gpd\n",
    "from matplotlib import rcParams\n",
    "\n",
    "# 设置字体为Arial\n",
    "rcParams['font.family'] = 'Arial'\n",
    "\n",
    "# 选择一个“Residential”节点及其对应的多边形\n",
    "sample_residential = residential_gdf.iloc[res_id]\n",
    "source_node = sample_residential['node']\n",
    "\n",
    "# 获取源多边形\n",
    "source_polygon = sample_residential['geometry']\n",
    "\n",
    "# 创建绘图对象\n",
    "fig, ax = plt.subplots(figsize=(12, 12))\n",
    "\n",
    "# 绘制道路网络\n",
    "ox.plot_graph(G, ax=ax, node_size=0, edge_color='lightgrey', show=False, close=False)\n",
    "\n",
    "# 绘制源多边形\n",
    "gpd.GeoSeries([source_polygon]).plot(ax=ax, color='#fcecb6', alpha=1, edgecolor='white')\n",
    "\n",
    "# 绘制源节点\n",
    "ax.scatter(sample_residential['node_x'], sample_residential['node_y'], color='red', s=50, label='Residential Node')\n",
    "\n",
    "# 标记是否已经绘制过某个目标类型的区域和节点，以避免重复图例\n",
    "area_labels_drawn = set()\n",
    "node_labels_drawn = set()\n",
    "path_labels_drawn = set()\n",
    "\n",
    "# 对于每个目标类型，绘制到最近目标节点的最短路径，并绘制目标多边形\n",
    "for idx, (t_type, target_gdf) in enumerate(target_gdfs.items()):\n",
    "    target_nodes = target_gdf['node'].tolist()\n",
    "    \n",
    "    # 计算到所有目标节点的最短路径长度\n",
    "    lengths = nx.single_source_dijkstra_path_length(G, source_node, weight='travel_time')\n",
    "    \n",
    "    # 找到最近的目标节点\n",
    "    target_nodes_in_lengths = [node for node in target_nodes if node in lengths]\n",
    "    if not target_nodes_in_lengths:\n",
    "        continue\n",
    "    nearest_target_node = min(target_nodes_in_lengths, key=lambda node: lengths[node])\n",
    "    \n",
    "    # 获取目标节点所在的多边形\n",
    "    target_row = target_gdf[target_gdf['node'] == nearest_target_node].iloc[0]\n",
    "    nearest_target_polygon = target_row['geometry']\n",
    "    \n",
    "    # 获取目标节点坐标\n",
    "    target_x = G.nodes[nearest_target_node]['x']\n",
    "    target_y = G.nodes[nearest_target_node]['y']\n",
    "    \n",
    "    # 绘制目标多边形\n",
    "    area_label = None\n",
    "    if t_type not in area_labels_drawn:\n",
    "        area_label = f'{t_type} Area'\n",
    "        area_labels_drawn.add(t_type)\n",
    "    gpd.GeoSeries([nearest_target_polygon]).plot(ax=ax, color=color_map[t_type], alpha=1, edgecolor='white', label=area_label)\n",
    "    \n",
    "    # 绘制目标节点\n",
    "    node_label = None\n",
    "    if t_type not in node_labels_drawn:\n",
    "        node_label = f'{t_type} Node'\n",
    "        node_labels_drawn.add(t_type)\n",
    "    ax.scatter(target_x, target_y, color=\"blue\", s=40, label=node_label)\n",
    "    \n",
    "    # 计算最短路径\n",
    "    path = nx.shortest_path(G, source_node, nearest_target_node, weight='travel_time')\n",
    "    \n",
    "    # 获取路径的坐标\n",
    "    path_coords = [(G.nodes[node]['x'], G.nodes[node]['y']) for node in path]\n",
    "    xs, ys = zip(*path_coords)\n",
    "    path_label = None\n",
    "    if t_type not in path_labels_drawn:\n",
    "        path_label = f'Path to {t_type}'\n",
    "        path_labels_drawn.add(t_type)\n",
    "    ax.plot(xs, ys, color=\"black\", linewidth=1, alpha=1, label=path_label)\n",
    "    \n",
    "# 获取当前图例的句柄和标签\n",
    "handles, labels = ax.get_legend_handles_labels()\n",
    "\n",
    "# 去重，保持顺序不变\n",
    "unique_labels = dict(zip(labels, handles))\n",
    "\n",
    "# 添加去重后的图例\n",
    "# ax.legend(unique_labels.values(), unique_labels.keys())\n",
    "\n",
    "# 添加标题并关闭坐标轴\n",
    "ax.set_title('Shortest Paths from a Residential Node to Target Nodes with Polygons')\n",
    "ax.axis('off')\n",
    "\n",
    "# 保存为PDF格式，确保矢量化输出\n",
    "plt.savefig('shortest_paths_to_targets.pdf', format='pdf', bbox_inches='tight')\n",
    "\n",
    "# 显示图形\n",
    "plt.show()\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 初始化结果字典\n",
    "average_travel_times = {}\n",
    "\n",
    "# 对于每个目标类型，计算平均旅行时间\n",
    "for t_type, target_gdf in target_gdfs.items():\n",
    "    print(f\"Calculating average travel time to {t_type}...\")\n",
    "    \n",
    "    # 获取目标节点集合\n",
    "    target_nodes_set = set(target_gdf['node'].unique())\n",
    "    \n",
    "    # 定义函数计算平均旅行时间\n",
    "    def compute_average_travel_time(G, source_node, target_nodes_set):\n",
    "        # 使用Dijkstra算法计算从源节点到所有可达节点的最短路径长度\n",
    "        lengths = nx.single_source_dijkstra_path_length(G, source_node, weight='travel_time')\n",
    "    \n",
    "        # 提取到目标节点的旅行时间\n",
    "        travel_times = [lengths[node] for node in target_nodes_set if node in lengths]\n",
    "    \n",
    "        if travel_times:\n",
    "            average_time = sum(travel_times) / len(travel_times)\n",
    "            return average_time\n",
    "        else:\n",
    "            return None  # 没有找到路径\n",
    "    \n",
    "    # 对每个“Residential”节点计算平均旅行时间\n",
    "    residential_gdf[f'avg_time_to_{t_type}'] = residential_gdf['node'].apply(\n",
    "        lambda node: compute_average_travel_time(G, node, target_nodes_set)\n",
    "    )\n",
    "    \n",
    "    # 将旅行时间从秒转换为分钟\n",
    "    residential_gdf[f'avg_time_to_{t_type}_minutes'] = residential_gdf[f'avg_time_to_{t_type}'] / 60.0\n",
    "    \n",
    "    # 处理缺失值\n",
    "    residential_gdf[f'avg_time_to_{t_type}_minutes'] = residential_gdf[f'avg_time_to_{t_type}_minutes'].fillna(-1)\n",
    "    \n",
    "    # 保存结果\n",
    "    average_travel_times[t_type] = residential_gdf[f'avg_time_to_{t_type}_minutes']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 计算到所有目标类型的总体平均旅行时间\n",
    "def compute_overall_average(row):\n",
    "    times = []\n",
    "    for t_type in target_types:\n",
    "        column_name = f'avg_time_to_{t_type}_minutes'\n",
    "        time = row.get(column_name, None)\n",
    "        if time is not None and time != -1:\n",
    "            times.append(time)\n",
    "    if times:\n",
    "        return sum(times) / len(times)\n",
    "    else:\n",
    "        return None\n",
    "\n",
    "residential_gdf['overall_average_travel_time_minutes'] = residential_gdf.apply(compute_overall_average, axis=1)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "from matplotlib import rcParams\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "\n",
    "# 设置字体为Arial\n",
    "rcParams['font.family'] = 'Arial'\n",
    "\n",
    "# 准备颜色映射：使用Blues，截取5段色阶\n",
    "cmap_name = 'OrRd'\n",
    "num_colors = 5\n",
    "cmap = mpl.cm.get_cmap(cmap_name, num_colors)\n",
    "\n",
    "for idx, t_type in enumerate(target_types):\n",
    "    column_name = f'avg_time_to_{t_type}_minutes'\n",
    "    \n",
    "    if column_name in residential_gdf.columns:\n",
    "        data = residential_gdf[column_name].dropna()\n",
    "        \n",
    "        # 计算各分位点：0%, 20%, 40%, 60%, 80%, 100%\n",
    "        quantiles = [0, 15, 30, 45, 60, 75]\n",
    "        \n",
    "        # 创建图形和坐标轴\n",
    "        fig, ax = plt.subplots(figsize=(12, 12))\n",
    "        \n",
    "        # 使用BoundaryNorm建立分段颜色映射\n",
    "        norm = mpl.colors.BoundaryNorm(quantiles, ncolors=num_colors, clip=True)\n",
    "        \n",
    "        # 绘制地图，使用离散的分区颜色\n",
    "        residential_gdf.plot(column=column_name, ax=ax, cmap=cmap, norm=norm, edgecolor='white')\n",
    "        \n",
    "        # 设置标题和字体大小\n",
    "        ax.set_title(f'Average Travel Time to {t_type} (minutes)', fontsize=20)\n",
    "        \n",
    "        # 创建colorbar\n",
    "        sm = mpl.cm.ScalarMappable(cmap=cmap, norm=norm)\n",
    "        sm.set_array([])  # 没有直接使用数组\n",
    "        \n",
    "        cbar = fig.colorbar(sm, ax=ax, fraction=0.046, pad=0.04)\n",
    "        cbar.ax.tick_params(labelsize=20, length=10)\n",
    "        \n",
    "        # 为colorbar设置刻度\n",
    "        # 这里会在colorbar上显示 quantiles 对应的值\n",
    "        # 如果想显示区间百分比，而不是对应的原始数据值，可手动设置cbar的tick labels\n",
    "        cbar.set_ticks(quantiles)\n",
    "        # 下面设置百分比标签（可根据需要调整格式）\n",
    "        percent_labels = [f'{int(q)}' if q > 0 else '0' for q in quantiles]\n",
    "        cbar.set_ticklabels(percent_labels)\n",
    "        \n",
    "        # 隐藏坐标轴\n",
    "        ax.axis('off')\n",
    "        \n",
    "        # 保存为PDF矢量格式\n",
    "        plt.savefig(f'avg_travel_time_to_{t_type}.pdf', format='pdf', bbox_inches='tight')\n",
    "        \n",
    "        # 显示图形\n",
    "        plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "from matplotlib import rcParams\n",
    "import numpy as np\n",
    "\n",
    "# 设置字体为Arial\n",
    "rcParams['font.family'] = 'Arial'\n",
    "\n",
    "# 列名称\n",
    "column_name = 'overall_average_travel_time_minutes'\n",
    "\n",
    "# 提取数据列并计算分位数\n",
    "data = residential_gdf[column_name].dropna()\n",
    "quantiles = [0, 15, 30, 45, 60, 75]\n",
    "\n",
    "# 创建离散颜色映射：这里用OrRd，并分为5个分级颜色\n",
    "cmap_name = 'OrRd'\n",
    "num_colors = 5\n",
    "cmap = mpl.cm.get_cmap(cmap_name, num_colors)\n",
    "norm = mpl.colors.BoundaryNorm(quantiles, ncolors=num_colors, clip=True)\n",
    "\n",
    "# 创建绘图对象\n",
    "fig, ax = plt.subplots(figsize=(12, 12))\n",
    "\n",
    "# 绘制分区地图\n",
    "residential_gdf.plot(column=column_name, ax=ax, cmap=cmap, norm=norm, edgecolor='white')\n",
    "\n",
    "# 设置标题和字体大小\n",
    "ax.set_title('Overall Average Travel Time from Residential Areas (minutes)', fontsize=19)\n",
    "\n",
    "# 创建colorbar，水平放置\n",
    "sm = mpl.cm.ScalarMappable(cmap=cmap, norm=norm)\n",
    "sm.set_array([])  # 不直接使用数组\n",
    "\n",
    "cbar = fig.colorbar(sm, ax=ax, fraction=0.046, pad=0.04, orientation='horizontal')\n",
    "cbar.set_label('Time (min)', fontsize=19)\n",
    "cbar.ax.tick_params(labelsize=19, length=10)  # 设置刻度字体大小和刻度线长度\n",
    "\n",
    "# 设置colorbar的刻度为分位数，并将其标注为百分比\n",
    "cbar.set_ticks(quantiles)\n",
    "percent_labels = quantiles\n",
    "cbar.set_ticklabels(percent_labels)\n",
    "\n",
    "# 隐藏坐标轴\n",
    "ax.axis('off')\n",
    "\n",
    "# 保存为PDF矢量格式\n",
    "plt.savefig('overall_avg_travel_time_horizontal_cbar_discrete.pdf', format='pdf', bbox_inches='tight')\n",
    "\n",
    "# 显示图形\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "# 保存结果到新的Shapefile或Feather文件\n",
    "\n",
    "residential_gdf.to_feather(\"residential_average_travel_times.feather\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "vector",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
